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ABSTRACT 

Procedures   for  designing   piecewise-linear    (PWL)    switching   surfaces 
for  minimum  time  control   are   investigated.      The  approach   is   to   generate 
the  switching   function  using  the   instantaneous   values  of   the  states   as 
inputs   to   easily  realized   PWL  function  generators.      The   design  problem 
is   then  one  of   fitting   the  state   space   PWL   surface  to   points   defining 
the  optimal    surface. 

Several   performance    indices   related   to   minimum-time  are   defined 
and   the  resulting   parameter   surfaces    for  different    initial   condition 
distributions   are  studied.      However,    because  of    local  minima,    a   success- 
ful  search   procedure   could   not  be   found.      Therefore,    a  heuristic  method 
of    least    squares   fitting   is  used.      A  combination  gradient   and  general- 
ized Newton-Raphson   search  method   is  employed   to   obtain  values    for  the 
PWL  parameters.      Several    least   squares    fit  methods   are   applied   to   a 
second  order   problem  and   the  results    compared   using   response   time  per- 
formance to   a  uniform  grid  of   initial   conditions.      Very  close   to   optimal 
performance    is   achieved   using   a  delay  time  design  procedure. 
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1.   INTRODUCTION 

This  study  is  concerned  with  the  specific  problem  of  minimum-time 
control  of  a  linear,  instantaneous,  time  invariant  plant  described  by  n 
first  order  linear,  constant  coefficient  differential  equation 

(1.1)       x  =  Ax  +  Bu 
where 

x  is  an  nxl  state  vector, 
A  is  a  known  nxn  plant  matrix, 

u  is  the  scalar  plant  control  which  is  bounded, 
B  is  a  known  nxl  distribution  matrix  of  the  control  u. 
The  plant  is  restricted  to  having  a  transfer  function  which  contains  no 
zeros  and  only  real  and  negative  eigenvalues.   In  addition,  the  n  states 
must  be  measureable. 

The  object  of  this  study  is  to  develop  a  general  design  procedure 
for  an  easily  "realizable"  controller  such  that  the  desired  minimum  time 
control  of  the  plant  can  be  effected.   A  block  diagram  of  the  system  is 
shown  in  Figure  1. 


>  controller 


7T 


cr 


contactor 


>  plant 


^ 


Figure  1.   Block  Diagram  of  Plant  and  Controller 
The  differential  equation  of  the  plant  can  be  normalized  with  respect 
to  the  magnitude  of  the  control  u,  such  that  I  u  1 5:  1.   The  control  u  is 
determined  by  the  scalar  switching  function,  <T  ,  such  that  u  =  sgn  C~. 
The  above  system  can  represent  either  (a)  the  so-called  "regulator" 


problem  where  the   input   c    is    zero   and  the  plant  must    be   forced    in  minimum 
time  from  some   non-zero,    random,    initial   condition   state   to   the  state 
space  origin,   or   (b)    the  general   control   problem  where  the    input    function 
c    is   non-zero,   but   restricted  to   an   allowable  class   of   functions,  e.g. 
step  or   ramp   for  a   second  order   plant.      This    implies   that   the   input   can 
be   imbedded   and  treated   as   an   initial   condition  of   the    error   state   vec- 
tor or   more   easily  visualized   as   a   translation  of   the   state   space  origin 
to   the  desired  point    in  state   space  e.g.    a   position  or  velocity  con- 
trolled  system.      Again  the   purpose  of  the   controller    is  to   generate   the 
proper   scalar    cT   so   as   to    force  the   plant   from   its    initial   state  to   the 
desired   state    in  minimum  time. 

Since  there    is   a  preponderance  of    literature  on  this   subject   of 
minimum  time   control,    a  one   paragraph  condensation  of   previous   work  will 
be   attempted  by  giving  only  the  most   significant   contributions    and  those 
which  directly   affect   this   study.      Historically,    the  minimum-time   regu- 
lator problem    (for  bounded   control,    constant-coefficient,    second   and 
third  order  plants)    was   first    studied  by  McDonald   [lj,    followed   by 
Bogner  and   Kazda  [2],    and   Bushaw  [3].      Bellman,   Glicksbury  and  Gross    [4] 
established   the   following   important   results: 

1.  The  optimum  control    law  is  bang-bang,    i.e.    minimum  time   response 
is    achieved    if  the  magnitude  of  the   control   u   is    always   at    its  maximum 
value,     |u|    =   1. 

2.  For  any   initial   condition   in  the   state   space,    the   state  origin 
will   be  reached    in  minimum  time  with,    at   most,    (n-1)    reversals  of   the 
control   u.      This   result   applies   to   any   constant-coefficient  nfc  -order 
plant   provided   the   plant   eigenvalues   are   real    and  negative. 

The  maximum  principle  of  Pontryagin  [5]   also   shows    that   the  control   u 
should   always   be  at    its  maximum  value.      For  the   case  where  the   plant 
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transfer   function  has    zeros    in  the   left-half   plane,    Schmidt    [6]   has   shown 
that    (n-1)    reversals   are   still   required  but   the  optimum  solution   is   not 
bang-bang.      The  distinction  must   be  made   that   this   result    is   not  for  the 
criteria  of   returning  the   state   vector  x  to   zero    in  minimum  time  but 
rather   returning  the   plant  output   and    its   derivatives   to   zero    in  minimum 
time. 

Previous    results   of   optimal   controller   design  can  be  divided    into  two 
classes;    those  using   a   digital    computer   as    the   controller   with   the   com- 
puter  programmed    to    rapidly   compute    the   required     0"   as    a    function  of 
time  using   the  differential   equations   of    the    plant   and   the    specific 
initial    condition  given,     <r  =   fCx^jt),   [7,8];    or   those   using   some   elect- 
ronic device  that   computes    the   required    (T  as   a    function  of   the    instan- 
taneous  states,    which    implies   that   the   states   must    be   available   for   feed- 
back.     The    first    approach   has   more   often  been   used,    however,    it    requires 
an  on-line  digital   computer.      The   second   approach    is   more   feasible  but 
has    not    been  widely   used   because   of    the   complexity   of   the  optimal    switch- 
ing  function.      This   function  cannot,    in  general,    be  expressed    in   closed 
form   for   third   or   higher  order   plants.  [9]        It    should   be   noted   that    this 
optimal    switching   function,    for   higher   than   second   order   plants,    defines 
a   switching   surface    in   state   space  where   this   curved   surface   partitions 
the   space    into  a   u   =    +1   region   and   a  u  =   -1   region.      Because   this 
switching  function  has   such  a  complicated   dependence  on   the    states,    the 
emphasis    in  controller   design  has   been  on  quasi-optimal    controllers 
which   are  more  easily   implemented   but  give    less    than  optimal   perfor- 
mance.   [10,11] 

As   shown   in  Figure   1,    this    study  will    use   the   second  method  of   com- 
puting the.  control   as    a   function  of   the  states    and  then   realizing  a 


controller  that    is  quasi-optimal.      The  method,    to   be  developed   in  Section 
2,    uses   piecewise-linear   functions  of   the   state  variables   to   generate   the 
quasi-optimal   switching  function.      Section    3  gives   a  phase  plane   plot 
analysis  of   the   effect  of   a  PWL  switching   line.      The  difficulties   en- 
countered  using  an   actual    response   time  performance   index   are  discussed 
in  Section  4,    compared  to   the   non-response   time   indices   discussed    in 
Section  5.      Section   6  gives   the  multi-dimensional   cost   surface   iteration 
search  procedure  which    is   used.      The  concluding   Sections   apply  the  devel- 
oped methods   to   2nd   and    3rd  order   system  examples. 
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2.       STATEMENT   OF   THE   PROBLEM 

The  minimum-time  optimal    switching   surface  which  divides    the   state 
space   into    u   =   +1   and   u   =   -1   regions    can   be  defined   by  the  optimal 
switching    function 

(2.1)  cr*   =    f*(x)   =   0 

where   the   ,lr  denotes  optimal.      Since   this    function  cannot,    in  general, 
be  expressed    in  closed   form  for   an   n        order   system,    an   alternative    is 


i  un- 


required   to   define    the    surface.      Smith    [l l]    has    used    a   reverse   time   s: 
ulation  of    the   differential    equations   of    the   plant   wherein   the   simulated 
plant    is    started    at    the  origin,    numerical    integration    is    performed    with 
time-running   backward,    and   the    sign  of    the   control    u    is    switched    n-2 
times.      The   result    is    a   set    of   points   describing   one-half   the   optimal 
switching   surface   i.e.    the   locus   of   one-half  of  all   trajectories   which 
takes    the   plant   to   the   state    space   origin    in    (n-2)    or    fewer   switching 
of   the  control    u.      The  one-half   comes    about   since   there   are   two   terminal 
trajectories    that   take   the    plant    to    the   origin   without    a   change   of    u, 
one   with   u   =    +1    and   the   other   with   u   =    -1.       It    can   be    shown   that    since 
the    plant    is    linear    and  the   control    u    is    symmetrical,    the   optimal    switch- 
ing  surface    is    an  odd    function  of   the   state  variables    and   therefore  only 
one-half  of   the   surface   need   be   generated.      It    is   assumed   that   the  ex- 
pect i-d   maximum   values    of    tin'    stales    of    the    plant    art.-    known    a    priori    such 
that    a  boundary   for   each    state  may   be   assigned.      Hence,    the   optimal    sur- 
face   is   bounded   and   the   number  of   points   describing   the  surface   is   finite 

In  addition  to   the  restrictions    placed  on  the  plant    in   Section   1,    a 
further  requirement    is   that   the  optimal   switching   surface  must   be   single 
valued    in  at    least  one   state.      Smith  fi.ll  has    shown,    using   Jordan  canon- 
ical  form,    that   for   a   linear   plant  with  real   eigenvalues,    the   time- 
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optimal  surface  is  single  valued  in  the  "uncoupled"  states, e.g.  described 
by  a  differential  equation  that  is  independent  of  the  other  states.   For 
this  study,  let  x-  be  this  "uncoupled"  state  variable  hereafter  called 
the  independent  state  variable.   The  optimal  switching  surface  can  now 
be  defined  by 

(2.2)      <T*  =  f*(x-j)p  +  (x.)p  =  0   for  j  =  1,2, ...n 

J  *   i 

p  =  1,2, ...P 

where 

P  =  total    number  of   points   defining  one-half   the  optimal 
switching  surface 

p   =    index  of  P  points 

x  .   =  dependent   states 
J 

x.    =   independent   state 
The  object    in  this   study   is   to    somehow  fit   the  optimal   switching 
surface    (points)    given  by  (2.2)   with  a   piece-wise   linear    (hereafter  PWL) 
surface  described   by   an  array  of  weights   W  as    shown   in  Figure  2. 


"j^V  Wj,6 


XWJIT 


_yx     for    j  =  l,2.  .  .n 


Figure   2.      General    PWL  Switching  Function 
Although  the  PWL  equations  may  be  written   in   several   different   forms, 
the  break  point   to   break   point   form  was   chosen  because  of   the   ease  of 


12 


digital    implementation.      Letting)^  be  the   number  of    segments   per  de- 
pendent   state  variable,    and  K  the  number  of  weights    per  dependent   state 
variable,    then  K  =   2  7\  +   2   since  each   additional    segment   beyond   the  first 
is  defined  by  two   additional   weights.      The  total    number  of   segments    is 
y-^(n-l)    and   the   total    number   of   weights    is    K(n-l).      The   set   of   weights 
can  be   arranged    as    an  array  of   numbers   with  no   algebraic    properties. 

(2.3)        |Kl)  (Wl,2)  Wl,3    •     •     •    W1,K-1         (W1,kTJ 

J(W2,l)  (W2,2)  W2,3    •     •     •    W2,K-1         (W2,K)/for    J=1.2,...n 

W  =  <  I  j    /    i 

wj,k  f  k=l,2,...K 

K-l,l)        K-1,2)       *n-l,3  W„-lfK-i    (*hJ  ^Yl+2 

With   the   PWL  surface   so    formulated,    W-    i_  with   k      _      refers    to    x.    value 

j , k        even 

of  the  break  point  with  W-  ,  •,  being  the  corresponding  x-  value  of  the 

J  »K~  1  J 

break  point,  where  &-  refers  to  the  PWL  function  generator  output  for 
x.  input.   Certain  of  the  weights  indicated  by  parenthesis  are  assumed 
to  be  held  fixed  by  the  formulation  of  the  problem,  e.g.  f W .  ,  and 
(w-  o)  are  determined  by  the  desired  origin  neighborhood  to  be  discussed 
in  Section  k,    and  (w-  K)  by  the  x.  boundary.   Therefore,  the  number  of 
variable  weights  p  equals  the  total  weights  minus  the  fixed  weights, 
P  =  (n-l)(2)^  -1).   Thus,  the  design  procedure  must  provide  values  for 
O  variable  parameters  such  that  the  PWL  surface  is  by  some  criterion 
fitted  to  the  optimal  surface. 

Although  this  study  is  aimed  at  realization  of  the  controller,  the 
emphasis  will  be  on  the  required  design  procedures  using  digital  compu- 
ter simulation.   However,  the  actual  PWL  switching  functions  can  be 
easily  realized  using  resistors,  diodes,  and  betteries  as  shown  in  Fig- 
ure 3.   The  W.  ,       values  are  the  battery  voltages  and  the  W-  ,   ,, 
j,k  even  :  6  J,k  odd 
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values  are  determined  by  the  various  resistor  combinations.  For  example 
in  Figure  3,  if  x-  is  between  W-  ,,  and  W-  fi,  diodes  (a)  and  (b)  are  con- 
ducting  and  current  is  flowing  through  RQ  and  R^;  thus 


(2.4) 


*J  -   "Rf 


\*o      *i/  J 


W-  2    Wi  4 
fl  Ro     Rl 


Figure  3a  shows  a  PWL  function  generator  for  >)=  4  and  x,-  2:  0.   An  equi- 
valent network  with  reversed  diodes  and  batteries  is  required  for  x.^  0. 


R3 

r-^mu> — fl- 


WJ,8 


*2 


W 
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Rl  b     wp^ 

— wwma/v — M — f 


R0 

— vAAAAAAA- 
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^-h5 — f 


j,2         wj,4      wj,6        wj,8 


(a)    PWL  Function  Generator   for  x.>0      (b)    Input   x.   vs  output   x-    for  x.>0 
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(c)   Controller,    3rd-0rder  Plant 
Figure    3:      Realization  of   PWL  Switching   Functions 
The  PWL  switching   function,   where  the  symbol  A    denotes   PWL,    cannot 
be  expressed  explicitly   in  terms  of  W  but   rather   in   terms   of   the  PWL 
function  generator  output   x-: 
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n 

£  ~  r~    ,  X-  * 


(2.5)  <T    =   0    =    f(x,W)    =    x1-    -   y>      x 

j  =  l       J 


or 


n 

(2.6)  xi   =  2[   *j 

j  =  l 

Thus,    the  controller,    as    shown    in   Figure    3c,    can  be   simply  realized. 

The   basic   question  which   remains   and   which    is    the   heart    of    this 
paper   is   how  to   fit   the  PWL  surface   with    p    variables    to   the   P  points 
of   the  optimal    surface   so    as    to    achieve   minimum-time   trajectories.      Two 
approaches   are   possible;    (1)    a   brute    force,    straightforward   technique  of 
setting  the   parameters   so    as    to   minimize   the   actual    response   times   to   a 
typical   set  of    initial   conditions   as   discussed    in   Section  4;    or    (2) 
defining   some  other  non-response   time   performance    index   as   discussed    in 
Section   5. 

In  order   to    investigate   the   problem   as    stated,    some   type   of    compu- 
ter   simulation    is    required.      A  digital   model   has    been   used   throughout 
this   study,    rather  than  an   analog  or   hybrid   model,    simply  because  of 
the   number  of   calculations    required.      Appendix  A  contains    the   digital 
model   and  the   computer   programs. 
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3.       EFFECT  OF   A  PWL  SWITCHING  SURFACE 

In  order  to  better  understand   the  non-linearities    produced  by  PWL 

switching,    a  qualitative   study  will   be  made  of   a  second-order,    double 

2 
integrator    (1/s    )    plant.      The  state  equations   of    this   plant   can  be  ex- 
pressed  as    a   function  of   the   state  variables   x,    and  x~  ,    and   control    u. 
The  resulting  two   first-order  differential   equations   are: 


(3.1) 


xl  "  x2 


x„   =   u,  where  for  minimum  time  control  u  =  *    1 
These  equations  can  be  solved  using  conventional  methods  in  terms  of  a 
specific  initial  state,  x  (tQ), 

for  u  =  -1    xx(tf)  =  x1(tQ)  +  tx2(tQ)  -  %t2 


(3.2)  x2(tf)  =  x2(tQ)  -  t 


and   for   u  =   +1  x, (tf )    =   xn (t    )    +  tx9(t    )    +  %t2 


■lvuf        Aivuo 


'2vuo- 


x0(tf)    =  x0(t_)    +  t 


■2^f 


2v,-o 


For  different    initial    states,    the   possible  trajectories    in  the  two- 
dimensional   state   space   are  two   families  of  parabolas   as   shown  in 
Figure  4a  and  4b  and  given  by 


for  u  =   -1  x.    =   -%x0      +  C 


(3.3) 


and   for   u  =   +1 


1  -2       •    "1 

:2   =   +%x22   +  C2 


>xl 


>xl 


(a)   Trajectories   for  u=+l        (b)    Trajectories   for   u=-l         (c)    Minimum-time 

Trajectories 

2 
Figure  4.      1/s      Plant    State   Space  Trajectories 
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Using  Pontryagin's  maximum  principle,    the  well   known  minimum-time  trajec- 
tories  to   force  the   states   from  some  random   state  to   the  origin  are  the 
partial    parabolas   as   shown    in   Figure  4c.      The  optimal   switching   surface, 
or   line    in   this   case,    is   the   trajectory  passing   through  the  origin  with 
n-2    =  0   switchings   which    is    (3.3)   with  C,    =  Co   =  0. 

If   this   optimal    switching   line    is   now  replaced   by   a   PWL  switching 
line,    what    does    a   typical    trajectory    look    like?     The   answer    is    that    at 
each   intersection  of   a   trajectory  with  the    PWL   line,    the   possibility 
exists   that   relay   chatter  or   after-end   point   action    (hereafter   called 
bumping)will  occur.      First,    consider  the  one   segment    PWL   line    (OC)    in 
Figure   5a  where   the   segment  divides   the  two    families   of   parabolas.      This 
is    the    familar    linear   switching   problem   as   obtained   with   tachometer 
feedback. 


>*! 


(a)  One  segment  (b)  Several  segments 

Figure  5.   1/s2  Plant  PWL  Switching  Trajectories 
It  can  be  seen  from  Figure  5a  that  for  each  possible  segment  whose  slope 
crosses  (at  point  B)  the  optimal  line,  some  member  of  the  u  =  +1  family 
will  be  tangent  to  the  segment  (at  point  A).   Thus,  any  u  =  -1  trajectory 
which  intersects  the  PWL  segment: 

(1)  Between  points  0  and  A  will  bump  down  to  point  0. 
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(2)  Between  points   A  and  B  will  move  with  u  =   +1   to   some   place 
between  points   0   and  A  and  then  bump. 

(3)  Between  B  and   C  will  move  with  u  =   +1   past   the  origin. 

(4)  Exactly   at   point   B  will  move  with  u  =   +1   to   the  origin    (optimal) 
For   the  general   case  with   several    segments,    Figure   5b,    similar  bumping 
action  occurs.      Any  u  =   -1   trajectory,   which  intersects   a  PWL  segment 

at    a  point   closer  to   the  origin  than  the  tangency   point    (A,   C  or    E)  , 
will   bump  until  the   slope   is   changed  by   the   next   segment    (points   B   and 
D) .      For  example,    a  u  =   -1   trajectory   intersecting  at   point  E,   will   bump 
down   to   point  D,   move  on  a  u  =   +1  path  and   either   pass   point   B  or   inter- 
sect  the  BD  segment  between  C   and   B  and  bump  down  to   B,    and  again  move 
on  a  u  =   +1   path   past  the  origin.      Therefore,   for  a   large  number  of   ran- 
dom  initial   conditions,    approximately  one-half  of   the  trajectories  will 
bump  upon  the  first   intersection  with  the   PWL  switching   line.      The  re- 
sults  are  easily   extended   to   higher  order   systems   where  the   PWL   surface 
again  crosses   the  optimal   surface  at    least  once  per   segment. 

Although  this  result    is    in  itself  trivial,    it   does   point  out  the 
PWL   switching  non-linearities   which  make   impossible  any  mathematical 
optimization  of   the   general    ntn-order,    PWL,   true  minimum  response-time 
problem.      Furthermore,    since   this   bumping   is   so   frequent  with  PWL 
switching    (without  other  modifications),    care  must   be  taken   in  modeling 
this   system  to    insure  that   the  bump  time  of  the  model   accurately   repre- 
sents  that  of   the  actual   relay  or   contactor.   For   example,    the  bump  time 
of   a  digital   computer  design  model    is   dependent  on  the   integration  step 
size;   or  for  an  ideal  relay  simulation,    somehow  recognizing  when  bump 
is  occurring   and   changing  to   an   (n-1)    set  of  non-linear  differential 
equations   dependent  on  the   segment   slope.      The   equation  of   the   latter 
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case   for  the   several   segment    1/s     plant   can  be   solved   in  terms  of   the 
bump   time,    tu: 

(3-4)  tb    =    1/b    log    I"*    *lX^l 

D  La+b    Xl(t0) 


where 


t,     =   ideal    relay  bump   time  to  move   from  x,  (t    )    t< 
x^Ctf)    along   a  segment 


b   =  the  slope  of    the  segment 
a   =   the   x2    intercept   of   the   segment 
For   higher  order   systems,    however,    the   equations   become   very   complex. 

It    is   realized   that    bump  or   continuous   rapid   changing   of   the   con- 
trol   is    not   tolerable   in   many   physical    systems.      One   procedure   to   be 
discussed    in   Section    7   does   eliminate  this   undesirable   feature. 
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4.       RESPONSE   TIME   PERFORMANCE    INDICES 

In  order  to  quantitatively   evaluate   the  operation  of   any   particular 
controller,    a   scalar   performance  criterion  or   cost,    denoted   by   J,    is 
defined  whose  value   depends   on  the   particular   setting  of   the   variable 
parameters.      The   problem  of   fitting   P  optimum  points   with   a   PWL  function 
of     P   variables   then  reduces   to   finding   the  minimum  value  of   J,    i.e. 
searching  a   performance   surface    in  a     P  +   1  dimensional   space   for   the 
global  minimum  of   J.      To  obtain  minimum  time  performance,    minimize 


•f 

(4.1)  J  ■■  dt   =  tf   -  tQ 


This  is  the  time  to  go  from  some  initial  state  x(t  )  to  some  final  state 
x(tf),  which,  for  this  problem,  is  the  state  space  origin.   However, 
since  the  PWL  surface  crosses  the  optimal  surface  several  times,  depend- 
ing on  the  number  of  segments,  more  than  several  initial  conditions  are 
required  to  adequately  define  the  fitting  of  the  two  surfaces.   In 
effect,  the  fit  of  the  two  surfaces  is  optimized  only  for  the  initial 
conditions  used  in  evaluating  the  cost  function.   Thus,  the  general  con- 
tours and  specific  minimums  of  the  cost  surface  are  both  a  function  of 
the  variable  parameters  and  a  function  of  the  number  and  distribution 
of  the  initial  conditions.   Since  there  are  no  established  criteria  for 
such  an  index,  several  different  cost  functions  were  studied: 

1   — 

Jl  "  —  ^-  ^ 
m=l 

<*•»      J2  *  TT  51  <*m  "  O 

m=l 

M 

J3^5T(tm-02 

m=l 
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M 

J4 =  4- SI  (vt*  - 1) 

^  M     ^— ^        m      m 

m=  1 


J5    =   -TJ-SICVt*   -    !>: 


m=l 

where 

M  =   total   number  of   initial    conditions   considered 

m  =    index  of    initial   conditions    (m  =    1,2,....M) 

*  th 

t      =   optimal    response   time   for   the  m        initial    condition,    x^Ct    ) 

tm   =    PWL  response   time    for   the  m        initial    condition   and   a   given 
PWL  surface 

The   true   average   optimal    response   time    for   M    initial    conditions    can   then 

be  given   in   terms  of   J^   and   J„   of    (4.2) 

1        M 
<*-3)  Jo   =  —  SI         <  =  Jl    "  J2 

m=l 

Each  of   the    five  definitions    in   (4.2)    provide   a  different   weight    for 
the  distribution  of   the    initial   conditions,   and   hence,    a  different   glo- 
bal  minimum   for   a   given    set   of    initial   conditions.      For   example,    J, 
gives   the  most   weight   to    those    initial    conditions   which  have   the   largest 

A 

t,    i.e.    the  ones    starting   farthest    from  the  origin;    J2   and   J3   also  give 
more  weight   to  those   farthest    from  the  origin  but   have    the   advantage 
over  J}   of   at    least    indicating  how  close  to   the   optimal    is   the  PWL;    J^ 
and   Jc   give  equal   weight   to   each   initial   condition   and   would    seem  to   be 
the  most   desirable.      Also,    the   squaring   in  J3   and  J  5   gives   more  weight 
to   the   larger  deviations   between  t   and   t      and   results    in  a   smoother   cost 
function.      However,    any  choice  of  one  over   the  other  four    is   arbitrary. 
Thus,    it  will   be   left    as   a  designer's   choice   according   to   the   particular 
problem  being   solved, i.e.    the  expected  distribution  of   the   actual    initial 
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conditions.   Having  once  obtained  t  and  t*  for  a  particular  initial 
condition,  it  is  an  easy  matter  to  compute  all  five  of  the  cost  func- 
tions and  hence,  all  five  will  be  compared  in  this  paper.   Appendix  A 
contains  subroutine  "COSTS"  which  computes  the  five  cost  functions  given 
by  (4.2)  using  subroutines  "OPTRT"  and  "PWLRT"  which  compute  t*  and  t 
respectively  for  a  given  initial  condition. 

As  an  introductory  study  of  the  minimum-time  cost  surfaces,  the 
1/s   plant  described  in  Section  3  was  chosen  using  a  one  segment  ( P  =  1) 
PWL  surface.   The  points  (P  =  14)  used  to  define  the  optimal  switching 
surface  were  generated  at  .1  second  intervals  by  integrating  backwards 
in  time  from  the  origin.   Before  proceeding,  several  seemingly  arbitrary 
decisions  must  be  made.   First,  since  it  is  not  possible  in  an  engineer- 
ing sense  to  reduce  the  states  to  exactly  zero  using  the  PWL  scheme, 
some  origin  neighborhood  must  be  defined.   For  this  preliminary  study, 
a  tear  shaped  origin  neighborhood  was  used,  whose  boundaries  are  trans- 
verse trajectories,  equation  (3.3).   Secondly,  since  both  states  are 
uncoupled  and  single  valued,  either  state  could  be  chosen  as  the  in- 
dependent state,  therefore  let  x-  =  x„.   The  boundary  for  the  dependent 
state  was  then  selected  at  x,  =  1.   The  third  arbitrary  decision  in- 
volves choosing  the  number  and  distribution  in  state  space  of  the  initial 
conditions.   As  a  first  attempt,  it  was  considered  desirable  to  have 
initial  conditions  spaced  away  from  either  of  the  two  surfaces,  say  on 
the  Xi  axis,  and  also  to  have  each  initial  condition  affect  an  equal 
number  of  P  points.   Hence,  the  initial  conditions  were  generated  by 
choosing  every  third  point  defining  the  optimal  surface  and  then  inte- 
grating from  these  points  backward  in  time  until  intersection  with  the 
x^  axis  occurred.   The  resulting  four  initial  conditions  were  x-^  =  .09, 
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.36,    .81   and   1.44;    x2    =   0.      Figure   6   shows   the   fourth  quadrant  of   the 

state   plane  with  the   four    initial    conditions,   the  optimal   trajectories 

and    the  origin  neighborhood.      Also   shown   in  Figure   6  are  four   possible 

segment   slopes    with  weights:      W,         =  W  =  0;    W  =   -1.667,    -2.222, 

1 »  1     1>2        1»3 

-3.333,  -6.667  and  w^  ^  =  1,  where  the  slopes  pass  through  the  optimal 
switching  line  at  the  four  corresponding  trajectory  intersections. 

W, 

(w1  1)=(w1)2)=0 


Optimal 


-Xi=-X, 


Figure    6:      Second   Order  Example 
Since   the  cost   functions   for  this    simple  problem  are  only   two  di- 
mensional,   J  can  be   easily   plotted   as    a   function  of   the  one   variable 
weight,   W"i    3.      Using   the  digital  model    shown   in  Appendix  A  and   an   inte- 
gration step   size  of    .0001   seconds,    the  optimal    and   PWL    (with  variable 
W-i    o)    response  times    for  the    four   initial   conditions   were   computed. 
Figure   7    shows   the  resulting   cost   surfaces   for  the   five  cost    functions 
defined   by   (4.2).      The  relative  magnitudes   of  J  are   unimportant.      How- 
ever,   the   shapes  of   the   surfaces,    specifically   the  multi-inflection 
points    (local   minimums),    present   a  very  difficult   problem  for  cost 
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functions   of  greater   than   two   dimensions,  i.e.    more   than  one  variable 
parameter.      The   task  of  searching   such  a  multi-parameter  surface    in  the 
presence  of   noise  (or    local   minimums)    is    an  area  of   current    investiga- 
tion.     However,   other  than  a  completely   random   search   procedure,    no 
successful  method   has   yet  been  developed.     [12] 

Using  the   above   second-order  example  with  revisions,    the   following 
results   were  obtained.      First,    a  one-initial-condition   problem  has    a 
cost    function  with   more   than  one   minimum.      In   fact,    it   has    a  minimum  at 
each  weight   (segment   slope)   which  results    in  an  optimal   switching   for 
a   given   number  of    half    revolutions    around    the   origin.       Four    initial 
conditions   resulted    in  more    than  eight   minimums,    with  one  occurring 
at   each  of   the    four  weights    shown   in  Figure   6.      Fourteen  and   twenty-four 
initial    conditions    were   also   tried.      The  result   was    that   the  general 
shape  of   the   cost    surface  became  more    linear   as   the   number  of    initial 
conditions   were    increased   but   also,    the   number  of    local   minimums   greatly 
increased   and   were   more   closely    spaced,    although   not    individually   as 
deep.      Thus,    several   conclusions   are:    (a)    the   number   and  distribution 
of    initial   conditions   directly   affects    the   number   and   spacing  of    the 
cost    function    local   minimums,    (b)    a    cost    function   minimum   exists    at 
each  weight   which  makes   a   particular    initial   condition  trajectory  op- 
timum,   (c)   using   a    large   number  of    initial   conditions   does    not   eliminate 
the   local  minimums   but   rather   increases    them,    and   (d)    the  general   shape 
of    the  cost   surface  does    improve  as   the  number  of    initial   conditions    is 
increased.      Hence,    it   seems   feasible    that   an  iteration  procedure   could 
be  developed  with   a   large  enough   step  size  so    as    to    "step  over"  the 
closely  spaced  minimum  provided  enough   initial   conditions   were   used. 
However,   the  difficulty   for   such  a  multi-parameter   search   now  becomes 
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one  of  excessive  computer  time  in  evaluating  the  response  times  for  each 
perturbation  of  the  search. 

In  an  attempt  to  improve  the  shape  of  the  cost  surface,  several  of 
the  arbitrary  decisions  previously  made  were  modified  in  turn.   Several 
origin  neighborhood  shapes  were  studied  including  circle,  square  and 
diamond.   The  only  conclusion  reached  was  that  the  transverse  trajectory 
size  of  the  neighborhood  had  a  direct  effect  on  the  width  of  each  mini- 
mum, i.e.  a  larger  neighborhood  resulted  in  a  wider  minimum,  but  in  no 
way  changed  the  general  shape  of  the  cost  surface.   Secondly,  x  was 
made  the  independent  variable,  but  again,  no  significant  change  in  re- 
sults were  noted. 

Frederick  [l(T]  has  studied  the  same  general  problem  as  above  with 
the  restriction  that  the  variable  parameters  of  the  PWL  switching  func- 
tion be  limited  to  a  small  number  so  that  the  cost  surface  may  be  ran- 
domly searched.   Methods  are  given  for  selecting  the  initial  conditions 
such  that  a  searchable  cost  surface  is  obtained.   Also,  root-locus 
techniques  are  given  for  designing  linear  (one  segment  per  dependent 
state)  switching  functions  in  order  to  initiate  the  search  of  the  cost 
surface  and  to  then  determine  which  of  the  components  should  be  made 
PWL.   However,  the  methods  are  not  adaptable  to  a  multi-parameter  prob- 
lem because  of  the  dimensionality  and  local  minima  of  the  cost  surface. 

Thus,  since  the  problem  as  stated  in  Section  2  cannot  be  readily 
solved  by  the  direct  approach  of  minimizing  the  response  times  of  a 
random  set  of  initial  conditions,  a  heuristic  method  was  used. 


26 


5.      NON-RESPONSE  TIME    PERFORMANCE    INDICES 

The   problem  of   fitting   for  minimum  response  time  the   P  optimum 
surface   points  with   a  PWL   surface  of    O    variables   can  be   indirectly 
solved    if   the   performance  criteria  or   cost    index   is   defined   such  that 
the   resulting    p    dimensional    cost    surface   can   be    successfully   searched. 
The   previous   section   showed    the    local   minima  difficulties    which   arise 
when   the   cost    index    is    given   as    a    function  of   the    response   times    for   a 
random  distribution  of    initial    conditions.      This    section  will   define 
two  other   indices   which   are    independent  of   actual    response   times   since 
they   both   use   the   common    least    squares    fitting   technique.      Hopefully, 
these    indices   will    result    in   setting   the    p    parameters   of    the    PWL  sur- 
face   for  minimum   response   times.      However,    since   the   true  optimal   PWL 
parameters    cannot   be   easily   found    per  Section  4,    the   use  of    either  of 
these    least    squares    fit    indices    to   design   the   controller   must    be    justi- 
fied  by   the   actual    response   time   performance    (4.2)    of   the   resulting   con- 
troller. 

Smith    [13]  defined    a   least    squares    fit    (LSF)    index  which    in   the 
notation  of   this   paper,    is 

1' 

P: 

where  x^  is  given  by  (2.6). 

In  words,  this  is  the  sum  for  P  optimal  surface  points  of  the  square  of 
the  difference  between  the  optimal  surface  (point)  and  the  PWL  surface, 
measured  in  terms  of  the  independent  state, i.e.  the  square  of  the  x^ 
distance  between  the  two  surfaces  measured  from  each  P  point.   Substi- 
tuting in  (5.1)  for  x^  yields 


p  =  l  HJ 
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(5.2)      J6  = 


£  R  xh  -  (x^\ 


For  the  one-segment  example  of  Section  4,  the  cost  surface  obtained  with 
(5.2)  has  a  smooth  parabolic  shape.   One  other  advantage  of  using  the 
LSF  index  is  the  elimination  of  the  arbitrary  decisions  concerning  the 
origin  neighborhood  shape  and  the  number  and  distribution  of  the  initial 
conditions.   Hence,  the  solution  for  p  parameters  is  unique  if  given 
only  the  P  points,  the  boundaries,  the  independent  state,  and  the  number 
of  segments.   Appendix  A  contains  subroutine  "SUMSQ"  which  computes  the 
cost  given  by  (5.2). 

Using  the  index  of  (5.2),  Smith  [l3]  developed  a  non-iterative 
analytic  method  to  fit  the  P  point  optimal  surface  with  a  PWL  surface 
of  %(  P  +  1)  parameters.   With  reference  to  Figure  2,  Smiths'  method 
uses  fixed  break  points,  i.e.  W,-  ^  even  are  fixed  at  values  which  are 
evenly  spaced  between  the  origin  and  the  boundary.   The  remaining 
%(  P    +1)  weights  are  found  by  defining  %(  0   +   1)  linear-segment  non- 
linearities  which  are  functions  of  the  dependent  states.   The  weights 
are  then  calculated  using  matrix  algebra.   The  method,  as  developed  by 
Smith,  will  be  used  in  this  study  (see  computer  subroutine  "STARTW" 
in  Appendix  A)  to  obtain  the  starting  weights  for  the  parameter  surface 
search  procedure  to  be  developed  in  the  next  section. 

A  second  possible  LSF  index  is  that  of  using  "deviation  time", 
6  j-,  defined  as  the  transverse  trajectory  time  spacing  between  the  opti- 
mal and  the  PWL  surfaces.   The  cost  function  is  given  by: 

P 
(5.3)      J7  =  X  <£?)- 
p=l     '  F 
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This  deviation  time  can  be  found  by  using  the  appropriate  sign  of  control 
u,  and  integrating  with  forward  or  backward  time  from  each  P  point  until 
intersection  with  the  PWL  surface.   If  the  cost  function  (5.3)  is  to  be 
evaluated  a  large  number  of  times,  such  as  in  an  iteration  search  method, 
it  would  be  easier,  rather  than  integrating  each  time,  to  define  a  grid 
of  points  on  each  side  of  the  optimal  surface  where  the  points  are  spaced 
a  given  number  of  time  increments  apart.   Figure  8,  using  the  1/s^  plant, 
shows  such  a  grid  in  which  there  are  two  points  on  each  side  of  each  op- 
timal surface  point. 


Optimal 

Figure   8.      Deviation  Time 
The   deviation  time   can  then  be   rapidly   approximated   by   counting   the 
number  of    increments    between   the  optimal    surface    point    and    the    PWL 
surface   and    interpolating   for  the   last    partial    increment.      The  method 
of   using  deviation  time   is  easily   extended   to   higher-order   systems   and 
would,    at   first  glance,    appear  to  give   a  better   LSF   in  terms   of   actual 
response  time.      The  cost   surface  obtained   with   (5.3),    for   the  one   seg- 
ment  example  of    Section  4,    has   a   skewed   parabolic   shape.      Appendix  A 
contains   subroutine   "OPTGRID"  which  generates    a  grid  of   points   on  each 
side  of   the  optimal    surface   and   subroutine    "DEVT"  which  computes   the 
cost   given  by    (5.3). 
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6.       SEARCH   PROCEDURE   FOR  MULTI -DIMENSIONAL  COST   SURFACE 

The   task  of   searching   a  multi-dimensional   cost    surface  would   be  very 
difficult   for   a   surface  containing    local   minima.      However,   by  using  either 
of  the   least   squares   fit    indices   of   Section   5,    a  relatively   smooth  para- 
bolic  shaped   surface   is  obtained.      Hence,    a  generalized   Newton- Raphson 
iteration  procedure  has  been  developed  which  has   the  advantage  of  rapid 
convergence. 

The   cost    index  J    (Jg  or  Jy)    is   assumed  to  be  a   function  of     p    para- 
meters,  W, ,   W2 5   W3,    .    .    .    .   Wp      where  each   parameter    is  one  of   the  variable 
weights  W^   u.  given  by  equation    (2.3).      J   is   then  approximated   by   a  Taylor 
series   expansion   in  the  W  space  where  W   is   a    p    x   1   column  matrix  of   the 
variable  weights.      The  expansion   is   terminated  with  the   second  order  terms 
of   the  Taylor   series   and   is    performed   about   the   point   W°  which   is    a    p    x   1 
column  matrix  of   the  starting  values   of   the   p    weights.      For   the   first 
iteration,   W°    is   obtained  by   Smith's    analytical  method  of  fixed   break 
points    as   discussed   in  Section  5.      For  the  following   iterations,    v/3    is 
the  W*    value  of   the  preceeding   iteration.      The   expansion   is: 

(6.1)         J(W)    =   J  (WO)    +    &jL      (W-W°)    +  ^(W-W°)Tl|^i(W-W°) 


where 


frAT  i 
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=  Q-jO  =  p  x  ^     matrix  of  2nd  partials 
therefore 

(6.2)  J(W)  =  J(W°)  +  (dV^  (W-tf3)  +  ^(W-W°)T  Q^  (W-W°) 

The  procedure    is   to   first   evaluate   J(W°)    using    (5.2)   or    (5.3),    then  per- 
turb  W  enough  times   to   generate   enough  equations   so   that  £  and  Q   can   be 
determined. 

If   the   variables    W-    are   perturbed    by   a   small    increment  A    >    one   at    a 
time    (i   =    1,2,.    .    .    .  r>  ) ,    where   W^    is   the    itn  component   of  W,    then    letting: 

wi+   =  W9    +  A         and   Wi_   =  W<?   -   A 
each   pair   of    perturbations    produces,    from    (6.2),    two   equations 

(6.3)  J(Wi+)    =   J(W°)    +   ri    (Wi+    -  W?)    +   igqii    (W.+    -   W?)2 

=   J(vf)    +     A  rL  +  h  A2qiL 

and   J(Wi_)    =   J(W°)    -     A    ri    +   1^A2qii 

These  two    simultaneous    equations    can   then   be    solved   for    r-    and   q- • 
yielding 

(6.4)  rt   =  _-L.[j(Wi+)    -   J^.)] 

qU    =   — L-  f  J(Wi  +  )    +    J(Wi_)    -    2J(W°) 
Thus    if   J(W°)    is    known,    and    if   each   variable    is    perturbed    by    + A    and 
-  A    and   the  corresponding   two   costs   evaluated,    the   r  vector   and   the 
diagonal   terms  of  Q   can  be  evaluated    from   (6.4).      This    is    so   far   a   total 
of   2  0    perturbations.      Next   the  q.  .   cross    partial   terms   must   be  evaluated. 
From  physical   considerations,    it    is  obvious  that  Q    is    a   symmetric  matrix, 
i.e.    q,-^   =  q-  •,    and  hence  only  one-half  of   the   cross    partial    terms   need 
be   calculated.      If  each   pair  of   two   variables    is    perturbed   by   + A  ,    then 
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substituting    in    (6.2)   yields: 

(6.5)  J(Wi+   +  Wj  +  )    =  J (W°)    +   rt    (Wi+   -   W°   )    +  r-(W.      -  W°   ) 

+  %qii(Wi  -  w°_)  +  %qjj(Wj  -  w°)  +  qij(wi  -  w?)(wj  -  w°o 

This  equation  can  be  solved  for  q-  •  since  the  r  terms  and  the  q^  terms 
were  previously  found  by  (6.4).   Rearranging  (6.5) 


(6.6) 


ri  +  rj 


qij-^JCg^V  -J(*>)]  --L 

To  find  all  the  q^,-  terms,  — i —  (  p-  1)  perturbations  are  required.   Hence, 
the  total  perturbations  required  to  find  £  and  Qis2p+-(p-l)  = 

JL   (P+  3). 

2    V 

The  new  weights   can   now  be   calculated   by  differentiating    (6.2)   with 
respect   to   W,    setting  the   resulting  equation  equal   to   zero,    and   solving 
for  the  new  weights,   W   . 


(6.7) 


Z-^-4    -   r   +  Q(W'    -  W°)   =  0 

\aw/      -     v-       - 


(6.8)  W*    =   W°    -   Q   lr_ 

For  the   scalar  case    (   P  =    1),    equation    (6.8)   reduces   to   the   familar 
Newton-Raphson   iterative  equation 

_dj 
(6.9) 


w'    =   W°   - 


dW 
d2j 

dW2 


Figure   9a   shows   the    scalar   case  Newton-Raphson  procedure  wherein  the 
object    is   to    find  the  value  of  W  for  which  J  =   f(W)    is    at   an  extreme 
point,    in  this   case   a  minimum.      The   necessary   and   sufficient    conditions 
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for  J  to  be  at  a  minimum  extreme  point  are  that =  0  and  >  0, 

dW  dw2 

The   Newton-Raphson  procedure    is    as    follows: 

1.  obtain   a   starting  value  of   the  variable   =  W° 

2.  approximate  — - —  with   a   straight    line  which    is   tangent   to 


dW  dW 

at    the  point    W° 

3.      find  the  value  of   W  =   W*    where  the   tangent    line    intersects 
dJ 


dW 


=   0 


dJ 


k.      repeat    steps    2    and    3   using   w'    as    the   new  W°   until   — —  =   0 


dW 


Equations  (6.8)  and  (6.9)  are  steps  2  and  3  for  the  generalized  p  variable 
and  the  scalar  (  p  =  1)  cases  respectively. 


>  W 


>   W 


>   W 


(a)       w  (b) 

Figure   9.      Scalar   Newton-Raphson   Iterative   Procedure 

To    insure   convergence  of   the   Newton-Raphson  method,    the   cost   surface 

d2T 
must   be  of   convex  shape,    e.g.  —  „.?0  for   the   scalar   case,   or   similarily, 

dW^ 
the   starting  weights   W°  must   be  within  this   convex  region  as   shown   in 

Figure   9b.      For  the  case  where  W°    is   near   the   permissible   Newton-Raphson 
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2 
boundaries,  e.g.i_i positive  but  near  zero,  a  slight  modification  of  the 
dW2 

procedure  is  required.   Equation  (6.8)  for  the  general  case  can  be  re- 
written 

(6.10)      w'  =  W°  -  K  Q_1r 
where  K  is  a  scaling  factor  (0£  K^  1).   Hence,  at  the  end  of  each  itera- 
tion, let  K  =  1,  find  w'  from  (6.10),  compute  J(w' )  and  compare  it  with 
J(W°).   If  J(W')  is  greater  than  J(W°),  let  K  be  less  than  one  and  con- 
tinue to  repeat  the  above  i.e.  a  single  variable  K  search  procedure  (see 
Section  9),  until  the  optimum  step  size  (within  1%)  is  found.   At  the  end 
of  the  iterative  process  (convergence),  each  element  of  r  will  be  zero  if 
J(w')  is  at  an  extreme  point,  and  Q  will  be  a  positive  definite  matrix  if 
J(w')  is  at  a  minimum  extreme  point.   Appendix  A  contains  subroutine 
"FITW"  which,  given  an  initial  w,    will  iterate  to  a  minimum  in  J  using 
(6.4) ,  (6.6)  and  (6.10). 
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7.   COMPARISON  OF  SEVERAL  PWL  DESIGN  PROCEDURES  FOR  1/s2  PLANT 

The  results  of  the  previous  sections  will  be  combined  and  used  in 
this  section  in  developing  and  evaluating  several  possible  PWL  design 
procedures  for  a  simple  second-order  plant  with  three  PWL  segments.   The 
general  design  method  will  be  to  fit  the  PU.  surface  to  the  P  points  of 
the  optimal  surface  using  one  of  the  LSF  indices  and  the  generalized 
Newton-Raphson  search  procedure.   Having  thus  designed  the  controller, 
its  performance  will  then  be  evaluated  by  comparing  it  to  other  possible 
controllers  including  the  true  optimum. 

In  order  to  measure  the  performance  of  any  particular  controller  de- 
sign, the  response  time  indices  given  by  (4.2)  will  be  computed  using  a 
fixed  set  of  initial  conditions.   Smith  [l3j  chose  the  P  points  of  the 
optimal  surface  as  the  set  of  initial  conditions  and  computed  only  the 
J,  index.   A  more  accurate  measure  of  performance  would  be  the  response 
times  of  a  complete  grid  of  initial  conditions  having  uniform  distribu- 
tion and  covering  the  entire  state  space  out  to  the  selected  boundaries. 
Because  of  symmetry,  only  one-half  of  state  space  need  be  used.   This 
method  would  better  test  the  entire  PWL  surface  and  would  more  closely 
approach  the  design  objective  of  a  completely  random  distribution  of 
initial  conditions.   Subroutine  "RECGRID"  of  Appendix  A  generates  this 
grid  of  points  in  the  state  half  space.   For  the  purposes  of  this  section, 
J-^  through  Jr  given  by  (4.2)  will  be  computed  for  both  initial  condition 
distributions,  in  order  to  note  the  different  quantitive  results  obtained. 
For  comparison  of  different  controller  design  methods,  the  1/s^  plant 
discussed  in  Section  3  was  chosen  with  boundaries  at  |x.|  =  1.125  and 
I X2 I  =  1.5.   The  optimal  switching  line  passes  through  opposite  corners 
of  this  boundary.   An  integration  step  size  of  0.01  seconds  was  used,  and 
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the  optimal  line  was  defined  by  15  points  spaced  0.1  seconds  apart.   Each 
design  method  was  performed  twice,  using  either  the  Jg  or  J   indices. 
To  evaluate  the  performance  of  each  design  method,  the  response  time  in- 
dices were  computed  for  both  the  15  initial  condition  points  defining 
the  optimal  surface  and  for  90  initial  condition  points  uniformly  distrib- 
uted in  the  state  half  space.   The  origin  neighborhood  used  was  again  of 
tear  shape,  having  cross  trajectory  boundaries,  equation  (3.3)  with 
c1  =  c2  =  0.01. 

The  several  methods  employed  with  their  corresponding  subroutine 
calling  arguments  are  as  follows: 

(1)  A  one  segment  PWL  surface,  (  p  =  1),  fitted  to  the  P  optimum 
points  with  the  generalized  Newton-Raphson  search  procedure  or  equivalent - 
ly,  Smith's  method  Ql3].   (CALL  OPTSUR,  STARTW,  COSTS,  RECGRID,  COSTS). 

(2)  A  three  segment  PWL  surface  with  fixed  break  points,  (P=  3), 
fitted  with  either  of  the  two,  for  this  case,  equivalent  methods.   (CALL 
OPTSUR,  STARTW,  COSTS,  RECGRID,  COSTS). 

(3)  A  three  segment  PWL  surface  with  variable  break  points, 

(  p  =  5),  fitted  with  the  Newton-Raphson  search.   (CALL  OPTSUR,  STARTW, 
FITW,  COSTS,  RECGRID,  COSTS). 

(4)  Same  as  (3)  above  except  that  the  5  PWL  parameters  were 
fitted  by  plotting,  i.e.  better  guess. 

(5)  Same  as  (3)  above  except  that  in  computing  the  response 
time  performance  the  feedback  loop  was  opened  such  that  only  one  switch 
was  allowed  per  state  space  quadrant.   (MODIFY  SUBROUTINE  PWLRT) . 

(6)  Same  as  (3)  above  except  that  the  P  optimal  points  were 
shifted  .06  seconds  backward  on  the  transverse  trajectory  before  being 
fitted;  and  in  computing  the  response  time  performance,  the  switching  of 
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the  relay  was  in  each  instance  delayed  .06  seconds.   (CALL  OPTSUR,  FITW, 
DELAY,  FTTW,  COSTS,  RECGRID,  COSTS,  ALSO  MODIFY  SUBROUTINE  PWLRT). 
Figure  10  is  a  table  of  the  design  weights  and  the  response  time  perfor- 
mance costs.   Method  1  was  included  in  order  to  show  the  advantages  of 
PWL  switching  over  standard  linear  switching.   Method  2  is  the  present 
published  "state  of  the  art"  [l3l.   Method  3  shows  the  degree  of  improve- 
ment obtained  by  using  variable  break  points  and  the  iterative  procedure. 
Method  U,   which  is  not  adaptable  to  higher  order  systems,  shows  that  the 
best  possible  PWL  fit  has  not  been  obtained  using  the  LSF  indices.   Method 
5  is  also  not  adaptable  to  higher  order  systems  but  it  does  indicate  the 
close  to  optimal  performance  that  can  be  obtained  when  the  relay  chatter 
(bump)  is  eliminated.   Method  6  yields  the  best  results,  again  by  elimi- 
nating the  bump.   Its  adaptability  to  higher-order  systems  and  realiza- 
tion will  be  discussed  in  the  next  section. 

A  study  of  the  data  in  Figure  10  would  suggest  the  following  comments 
and  conclusions: 

(1)  Jg  versus  J7 : 

In  method    1,   J,   gives   a  better   fit   than  J-,    for   the  optimum  sur- 
face   initial   conditions    (I.C.'s)  whereas   J-,    is   better   than  J,    for   the 
uniform  grid  of   I.C.'s.      In  methods  2    and    3,   Jg  gives    a   better   fit    for   both 
initial    condition  distributions.      In   the   remaining   methods,    there    is    no 
comparison  between  Jfi   and  J-j   because   of   the  differences   between  J,    through 
Jr.      Hence,    the  conclusion   is   that  Jfi   and   J7   give   roughly  the  same   re- 
sults.     However,    since  Jft    is  much   easier  to   compute    (digitally   faster), 
it    is   preferred   over   the  more  difficult   J-,. 

(2)  Comparing  J-^  through  J5: 

There  is  no  apparent  "goodness"  of  any  one  over  the  four  since 
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each  provides  similar  information  although  weighting  the  individual  res- 
ponse times  differently.   As  stated  in  Section  4,  it  will  be  left  as  a 
designer's  choice  depending  on  the  expected  distribution  of  actual  initial 
conditions-   As  a  personal  choice,  Jk  will  be  used  in  the  remainder  of 
this  study  because  it  gives  equal  weight  to  each  initial  condition  and 
inherently  includes  the  comparison  with  the  true  optimal  i.e.  J  =  1  says 
that  for  a  given  initial  condition,  the  expected  value  of  t  will  be  twice 
that  of  tm. 

(3)  Optimum  surface  initial  conditions  versus  the  uniform  grid  initial 
conditions:   From  Figure  10,  using  the  optimum  surface  initial  conditions, 
the  costs  for  method  (1)  with  one  segment  appear  better  than  method  (2) 
with  three  segments,  which  is  obviously  incorrect.   The  error  comes  about 
because  all  the  optimum  surface  points  (I.C.'s)  which  are  closer  to  the 
origin  than  the  first  PWL  segment  intersection,  are  already  beyond  the 

PWL  switching  and  hence  are  optimum.   On  this  basis,  it  is  concluded  that 
the  uniform  grid  of  initial  conditions  better  tests  the  entire  PWL  switch- 
ing surface. 

(4)  Comparison  of  the  five  methods  yields  the  following  conclusions: 

a.  PWL  switching  is  better  than  linear  switching. 

b.  Optimizing  the  break  points  is  better  than  arbitrary  break 
points. 

c.  The  LSF  indices  do  not  provide  the  best  response  time  PWL 
fit. 

d.  Eliminating  the  bump  or  relay  chatter  produces  very  near 
o  pt  ima 1  per f o  rmance . 
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1/j2  PUNT 
Beundariea:  jx^ml.1^5,  foe   1=1.5 

DESIGN 


MfiTHODb 

WEIGHTS 

Fit 
Index 

(Parenthesis    indicatea    fixed) 

W 
*l-3 

-1.9*5 

*14 

l< 

W16 

"'7 

W18 

1.    linear 
segments 

J6 

(1.125) 

J7 

-^•559 

(1.125) 

2.    fixed 

break    pta 

n-5 

^*5 

J  ,5 

-1.020 

(0.575) 

-1.204 

(0.750) 

-1.551 

(1.125) 

J7 

-1.052 

(0.575) 

-1.190 

(0.750) 

-1.516 

(1.125) 

'.    variable 
break   pta 

e*5 

J6 

-O.456 

o.ei2 

-0.956 

O.456 

-1.526 

(1.125) 

J7 

-0.514 

0.512 

-0.646 

0.550 

-1.545 

(1.125) 

4,    same   as 
5.    abova 

graph- 
ical 

-0.450 

0.065 

-0.660 

0.400 

-1.550 

(1.125) 

^.    same   aa 
%    above 
with   open 
loop 

J6 

-O.456 

0.612 

-0.956 

O.456 

-1.528 

(1.125) 

J7 

-0.514 

0.512 

-0.648 

0.550 

-1.5*5 

(1.125) 

6.    same    as 
5.    abeve 
with   delay 
time 

J6 

-0.257 

0.526 

-0.797 

0.587 

-1.404 

(1.125) 

'7 

-O.256 

0.055 

-o.75£ 

0.556 

-1.416 

(1.125) 

LEGEND: 


Non-Response  Time    Indices 
J6  = 


%\%m-]n   '-sw. 


L    j\i 


FIGURE  1C.  Result?  of  Second  Order  PWL  Switching 
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FIGURE  10.   Cantinuad 


METHODS  ||  RESPONSE  TIKE  PERFORMANCE 

OptimuB 
=15  wit 

•urfaca  initial  canditiane   ]        Dnifara  grid   initial  aanditlant 
h     J.=  0.614                                            =90  with  J.=2.009 

'i 

J2 

J5 

J4 

J5 

Jl 

h 

J5 

h 

J5 

1. 

J6 

0.962 

0.348 

0.657 

0.559 

0.150 

3.084 

1.074 

1.546 

0.647 

O.605 

J7 

1.223 

O.609 

1.239 

1.034 

0.400 

2.955 

0.945 

1.222 

0.555 

0.445 

2. 

J6 

1.217 

0.605 

0.802 

0.955 

0.472 

2.905 

0.896 

1.038 

0.492 

O.56O 

Jy                1.242 

0.629 

0.872 

0.999 

0.496 

2.915 

0.905 

1.086 

0.495 

0.576 

5. 

J6    [    0.895 

0.281 

0.109 

0.403 

0.257 

2.324 

0.315 

0.119 

0.186 

0.074 

J7 

0.926 

0.312 

0.123 

0.605 

0.598 

2.340 

0.550 

0.151 

O.196 

0.077 

4. 
5. 

6. 

grap- 
hical 

0.854 

0.240 

0.079 

0.547 

0.225 

2.298 

0.292 

0.104 

0.175 

0.066 

J6 

0.727 

0.113 

0.042 

0.246 

0.043 

2.095 

0.087 

0.062 

O.058 

0.041 

J7 

0.703 

0.089 

0.051 

0.548 

0.052 

2.105 

0.093 

0.045 

0.051 

0.052 

J6 

0.791 

0.112 

0.043 

0.298 

0.011 

2.071 

0.06£ 

0.056 

0.052 

0.015 

J7 

0.806 

0.140 

O.O36 

0.289 

0.039 

2.056 

0.054 

0.049 

0.048 

0.022 

LEGEND: 


Raspansa  Tinaa   Indies*: 


m  tn 


2 


J1~V    -^-   tn 


■r*  £  k  -  d 

33*1    L  J 


"■it  ft-] 

wfc[H 
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8.   GENERAL  DESIGN  METHOD 

The  previous  section  indicated  the  desirable  results  obtained  with 
a  PWL  switching  function  and  the  delay  time  method  of  eliminating  the 
bump.   The  generalization  of  this  method  is  explained  in  this  section  and 
a  possible  delay  time  switching  circuit  is  given.   The  general  delay  time 
design  procedure  for  a  plant  of  order  n  is  as  follows: 

a.  Generate  the  P  points  defining  the  optimal  switching  surface  by 
starting  at  the  origin,  integrating  the  state  equations  with  reverse 

time  and  switching  the  control  n-2  times.   This  requires  defining  a  bound- 
ary (maximum  expected  value)  for  each  of  the  states,  an  integration  step 
size  and  the  time  spacing  between  the  P  points.   Subroutine  "OPTSUR"  in 
Appendix  A  generates  the  P  points. 

b.  Fit  the  P  points  of  the  optimal  surface  with  a  PWL  surface  of  P 
variables  using  a  least  squares  fit  cost  index  (J^  given  by  equation 
5.2)  and  an  iterative  procedure  (Newton-Raphson ,  equation  6.10)  to  search 
the  cost  surface.   It  requires  choosing  an  independent  (single  valued) 
state  and  the  number  of  PWL  segments  per  dependent  state.   Subroutine 
"FITW"  in  Appendix  A  fits  by  least  squares  the  PWL  surface  to  the  optimal 
surface  points,  i.e.  finds  values  for  /i  variables. 

c.  Find  the  maximum  deviation  time,  (  &  +.)    ,  between  the  optimal 

0   t  max'  r 

surface  and  the  PWL  surface  by  integrating  from  each  P  point,  with  for- 
ward or  reverse  time,  along  the  transverse  trajectory  until  intersection 
with  the  PWL  surface. 

d.  Shift  each  of  the  optimal  surface  points  2(  £   t^max  un*-ts>  back- 
ward in  time  along  its  transverse  trajectory.   Subroutine  "DELAY"  in 
Appendix  A  computes  (  S  t^max  an<^  shifts  the  P  optimal  surface  points. 

e.  Refit  the  P  points  of  the  shifted  optimal  surface  with  a  p  var- 
iable PWL  surface  as  in  step  b.  but  without  changing  the  independent  state 
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or  the  number  of  segments. 

f.   Realize  the  controller  with  (n-1)  PWL  function  generators,  Figure 
3,  using  the  p  weights  of  step  e.   Also  use  a  contactor  or  relay  with  an 
inherent  delay  of  2(&   t)    ,  i.e.  switches  the  control  u  at  2( S  t^max 
time  units  after  the  scalar  switching  function  <J~  passes  through  zero. 

The  above  can  be  most  easily  visualized  in  the  two  dimensional  state 
space  by  the  1/s^  example  of  Section  3.   Shown  in  Figure  11  are  the  P 
points  of  the  optimal  line,  the  PWL  fit,  the  maximum  deviation  time 
(  $  .-)    and  the  shifted  P  points.   The  (  i  t^max  t^-me   units  on  each  side 
of  the  original  optimal  line  defines  a  swath  2($  t)max  time  units  wide 
within  which  most  trajectories  will  remain  after  the  first  switch.   How- 
ever, note  that  the  boundaries  of  the  swath  do  not  define  trajectories 
but  rather  in  moving  toward  the  origin  the  swath  tends  to  get  closer  to 
the  optimal  whereas  the  trajectories  tend  to  diverge.   Hence,  not  all 
trajectories  will  remain  in  the  swath  after  the  first  switch.   If  each 
P  point  of  the  optimal  line  is  shifted  backward  in  time,  2(  £  t)    units 
along  its  transverse  trajectory,  a  second  swath  is  obtained.   Now,  if 
the  switching  function  a"   goes  through  zero  in  the  shifted  swath  (at  the 
PWL  surface  which  has  been  fitted  to  the  shifted  P  points)  and  if  the 
actual  relay  switching  is  delayed  2(6t)    ,  then  the  actual  switching 
will  occur  in  the  original  swath.   Thus,  the  objective  of  eliminating 
the  bump  will  be  completely  achieved  for  all  trajectories  remaining  in 
the  swath  and  partially  achieved  for  those  leaving  the  swath. 

The  delay  time  method  requires  a  relay  that  switches  the  control  u 
a  fixed  time  after  the  switching  function  a~   goes  through  zero.   One 
awkward  but  workable  device  to  accomplish  this  is  shown  in  Figure  12 
using  three  relays  and  an  RC  changing  circuit.   Whenever  ff"  passes 
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Figure    11.       Delay   Time  Method 
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Figure  12.   Delay  Time  Switching  Circuit 
through  zero,  relay  1  switches  and  charges  capacitor  C  with  a  D.C.  volt- 
age t   V  which  in  turn  energizes  relay  2  and  opens  the  input  to  relay  3. 
This  relay  3  with  feedback  has  the  overall  characteristic  shown  in  Figure 
12  where,  in  this  case,  the  constant  K  is  just  large  enough  such  that  it 
does  not  switch  when  relay  1  does.   The  RC  discharge  time  and  the  drop 
out  current  of  relay  2  determines  the  delay  time  as  shown  in  Figure  12. 
Hence,  at  the  end  of  the  delay  time,  relay  2  closes  and  since  fl"  is  now 
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greater  than  K,  relay  3  switches  the  control  u.   The  above  time  delay 
switching  could  probably  be  more  elaborately  done  with  digital  circuitry 
using  a  Schmidt  trigger  arrangement.   However,  the  above  is  satisfactory 
for  indicating  the  realization  of  time  delay  switching. 

The  results  of  using  the  proposed  time  delay  method  are  indicated  in 
the  table  of  Section  7  for  the  second  order  plant.   For  all  intents  and 
purposes,  the  response  times  achieved  with  PWL  switching  and  delay  time 
are  within  some  origin  neighborhood  of  the  true  optimal.   The  next  section 
will  attempt  to  apply  the  proposed  method  to  a  third-order  example. 
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9.      THIRD  ORDER   EXAMPLE 

The  previous  section  defined  the  general  design  procedure  which 
will  be  applied  in  this  section  to  a  triple  integrator  (1/s^)  plant. 
The   state   equations    for   this    plant   are: 


(9.1) 


xl    "   x2 


x2    ~    x3 


Xo   =   u        where   u   =    t    1 


The   time   optimal    switching   surface    for   forcing   the    states    to    zero    from 
any    initial    condition   can  be   expressed    in   closed    form    (_14J    as: 


(9.2) 


3 
*  xl 

<r    ~~  °  =  xi  +  — r  + 


x23 


G  ^x,  +  £( — 2-  +ex9) 


3/2 


where 


€=   sgn 


X2    + 


x3    !x3l 


The   resulting   surface    is   singled   valued    in   all    three    states    and    is    sym- 
metrical   about   the   state   space   origin   since   the   control    u    is    symmetrical. 

The   design  objective   for   this   example   will   be   to    approximate   the 
optimal    surface  with   a    PWL   surface  of    four    segments    (  Y[  -   4)    per   depen- 


^  5l  , 


^  II 


^  6. 


dent  state,  within  the  region  bounded  by  x. 

From  Section  2,  the  number  of  variable  weights  is  P  =  (n-l)(2)1  -1)  =  14. 

The  other  rather  arbitrary  choices  are:   let  x,  be  the  independent  state, 

let  the  integration  step  size  be  0.1  seconds  and  let  the  time  spacing 

between  the  optimal  surface  points  be  1  second.   Using  Subroutine 

"OPTSUR"  in  Appendix  A,  the  optimal  half  surface  is  then  defined  by 

30  points. 

The  second  step  outlined  in  the  procedure  of  Section  8  is  to  now 
search  the  f   +  1  =  15  dimensional  space  for  the  minimum  LSF  cost  as 
given  by  equation  (5.2).   The  search  procedure  is  the  generalized  Newton- 
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Raphson   iterative  method  of   Section    6.      However,    because   some  of   the 
initial   values   W    of   the   14  variables    are  apparently  not   within  the 
allowable  Newton-Raphson   range    (see  Figure   9),   the  method  was   found 
not   to  converge.      Also,    because  there   exists   such  a  complicated   relation- 
ship between  the  variable   weights,    the  cost    surface  was   found   to   contain 
many   irregularities    i.e.    inflection   points,    saddle   points,    etc,   which 
hindered  the   search. 

Since   Smith's  method  £l3j    is   the  only  reasonable  starting   point 
for  W°,   a  different   starting  search   procedure   is   required.      Using  equa- 
tion  (6.4)    for  r    (  p  x   1  matrix  of   first   partials),    a  gradient  or 
steepest   descent  method   can   be  defined   as: 

(9.3)  w'    =   if    +   Kr 

where  K  is  a  constant  defining  the  step  size. 

The  difficulty  with  the  gradient  method  is  in  choosing  a  proper  K  since 
too  small  a  step  will  result  in  a  very  slow  descent  (many  iterations) 
and  too  large  a  value  could  result  in  an  increasing  cost  (an  ascent)  or 
an  overlapping  of  the  PWL  break  points.   One  possible  definition  for 
the  step  size  is: 

(9.4)  K  =   AJ 


'r_Tr 
where  A  J  =  J(W' )  -  J(W°)  =  the  change  in  cost  desired  which  is  defined 

a  priori  (before  the  step) . 

For  the  purpose  of  this  study,  each  computation  of  r  will  be  con- 
sidered as  one  gradient  iteration.   Then,  for  each  gradient  iteration, 
a  single  variable  search  procedure  was  developed  to  find  the  value  of 
K  which  minimized  J(w' ) .   The  reasons  being  that  J(w')  was  found  to  vary 
significantly  with  the  choice  of  K  and  secondly  each  gradient  iteration 
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requires    2f>  (28    in   the   third   order   example)    computations   of   the   cost 
versus   one   computation    for  one   value   of    K.       Each   gradient    iteration    is 
then   analogous   to   determining   the   direction  of    steepest    descent    r    from 
some   point   W°   and   then  moving    in   this  direction   until    no    further   descent 
is    possible  without   changing   directions.      The    single   variable   search 
procedure    is    as    follows: 

1.  Compute   the   starting   value   of    K  from    (9.4)   with     £j  equal    to 
one    per   cent   of   the   preceeding   J(w'). 

2.  Compute   w'    from    (9.3)    and   the   resulting   J(w'  )  . 

3.  a.    If   J(w')<    J(W°),    let    J(W°)    =   J(W'),    and    repeat    step   2. 

b.    or    if   J(W')   >  J(W°),    let    K   =    -.IK,    J(W°)    =   J(W),    and   repeat 
step   2. 

The   process    is    continued    for   three    reductions   of    K,    i.e.    K    increasing 
in    increments   of    1,    decreasing    in    increments   of    .1,    and    increasing    in 
increments    of    .01,    such   that    the    final    value   of    K    is   within   one   per 
cent   of    the   optimal    value. 

The    above  gradient   method   was    tried   on   the     p    =    14,    third— order 
problem,    but    again   because  of    the  many    irregularities   of   the    cost    sur- 
face,   the  method   did    not    converge   after   numerous    iterations,    although 
each    iteration  did   result    in   an    improvement.      A   final    modification  was 
that   of    a   hybrid    search   where   the    gradient   method   was    used    for   start- 
ing.     After   each  gradient    iteration,    the  q^    (self    second    pcirtials)    were 
computed,    and    if   all    q--    were    positive    (a   necessary  but    not    sufficient 
condition    for  Q   to   be    positive   definite),    then  the   search   automatically 
switched   to   the   Newton-Raphson  method   given   in   Section   6.      This  hybrid 
method   also    failed   because  the   test   of  Q   for  going   to   Newton-Raphson 
was    not   sufficient   to    insure   the   step  was    in  the  proper  direction   i.e. 
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Q  positive  definite. 

The  two  questions  of  whether  or  not  some  other  iterative  search 
procedure  can  locate  the  global  minimum  of  J  (fit  the  PWL  surface)   in 
a  reasonable  amount  of  computer  time  and  whether  or  not  the  delay  time 
method  of  eliminating  the  bump  is  adaptable  to  third  and  higher-order 
systems  will  be  the  subject  of  further  study. 
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10.   CONCLUSIONS  AND  ACKNOWLEDGEMENTS 

The  problem  of  designing  an  easily  realizable  minimum  time  (bang- 
bang)  controller  was  investigated.   The  approach  was  to  realize  the 
switching  function  by  using  the  instantaneous  values  of  the  states  as 
inputs  to  piecewise  linear  function  generators.   The  problem  then  re- 
duces to  finding  the  optimal  weights  for  the  piecewise  linear  functions 
i.e.  fitting  the  PWL  surface  to  the  optimal  surface. 

A  qualitative  study  of  the  effect  of  PWL  switching  showed  the  un- 
desirable but  very  prevalent  relay  chatter  (after  end  point  action) 
which  results.   Several  minimum  time  performance  indices  were  defined 
and  the  resulting  performance  surfaces  for  different  initial  condition 
distributions  were  studied.   The  conclusion  reached  was  that  because 
of  local  minima,  the  global  minimum  of  these  performance  surfaces  could 
not  be  readily  found.   Therefore,  a  heuristic  design  method  was  re- 
quired, namely,  least  squares  fitting.   The  search  of  the  least  squares 
performance  surfaces,  for  the  multi-dimensional,  second-order  case, 
converged  quite  rapidly  using  a  generalized  Newton-Raphson  iterative 
procedure. 

Several  different  design  methods  for  a  second-order  plant  were 
studied,  including:   linear  switching,  arbitrarily  fixed  PWL  break 
points  [13],  variable  break  points,  and  delay  time.   It  was  concluded 
that  PWL  switching  is  better  than  linear  switching,  that  variable  break 
points  is  better  than  fixed  break  points,  and  that  delay  time  is  very 
near  the  true  optimal  since  it  eliminates  the  relay  chatter  difficulty. 
This  comparison  of  methods  was  based  on  actual  response  time  performance 
to  a  uniform  grid  of  initial  conditions.   It  was  also  concluded  that  a 
uniform  grid  of  initial  conditions  in  the  state  space  better  tested  the 
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performance  of  a  PWL  switching  surface  than  did  an  optimum  surface  grid. 

The  variable  break  point  and  delay  time  methods  were  then  applied 
to  a  third  order  plant  with  fourteen  PWL  variables.   The  Newton-Raphson 
procedure  for  fitting  the  PWL  surface  would  not  converge  because  of  the 
starting  values  of  the  parameters.   The  use  of  a  gradient  (steepest 
descent)  method  was  also  tried  but  did  not  converge  after  numerous 
iterations  although  each  iteration  did  result  in  an  improvement. 

The  problems  for  further  study  are  numerous.   Among  the  most  impor- 
tant are  a  search  procedure  which  converges  in  a  reasonable  time  for 
third  and  higher-order  plants  and  a  method  (perhaps  delay  time)  which 
eliminates  the  relay  chatter  for  third  and  higher-order  plants.   There 
are  many  facets  of  this  eliminating  the  bump  problem  which  have  not 
been  fully  investigated  in  this  paper.   One  seemingly  easier  approach 
than  delay  time  would  be  that  of  switching  on  some  constant  tr  value 
as  shown  for  relay  #3  (with  feedback)  of  Figure  12.   For  the  second- 
order  case,  the  constant  <r*  loci  are  translations  of  the  <T  =  0  (PWL) 
line  in  the  independent  state  direction.   However,  for  a  third-order 
system,  the  information  for  the  second  switch  from  the  optimal  surface 
to  the  optimal  line  would  seem  to  be  lost.   Another  interesting  problem 
that  could  be  studied  is  the  relationship  between  the  desired  origin 
neighborhood  size  and  the  minimum  number  of  PWu  segments  required  to 
avoid  overshooting  the  neighborhood.   It  would  appear  that  the  width 
of  the  maximum  deviation  time  swath  around  the  origin  would  define  a 
desirable  origin  neighborhood  size  and  hence  indirectly  relate  the  num- 
ber of  segments  required  for  a  given  origin  neighborhood.   Other  refine- 
ments include:   using  linear  combinations  of  the  states  to  produce  cross 
product  terms  (curved  PWL  surfaces)  [l3j;  and  possibly,  some  other  sur- 
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face   fit   criteria.      A  close   to   optimal   general   problem  solution  would 
be  of   considerable    interest   since   it   would   help  bridge   the  gap   between 
optimal    theory   and   optimal    practice. 

The    author   would    like   to   express  his    appreciation   for   the   guidance, 
assistance   and   encouragement   of   Professor   D.    E.    Kirk   under  whose  direc- 
tion  this    study   was    performed.       In   addition,    thanks    are   due    to    the 
author's    fellow  graduate   students    for  many   helpful   discussions. 
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performance    is    achieved   using   a  delay   time   design   procedure. 
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